function C = covariance(varargin)
    C = zeros(nargin);
    for i=1:nargin
        for j=1:nargin
            tmp = (varargin{i}/PhysicalVal(1,0,varargin{i}.units))+(varargin{j}/PhysicalVal(1,0,varargin{j}.units));
            C(i,j) = (tmp.stderr^2 - varargin{i}.stderr^2 - varargin{j}.stderr^2)/2;
        end
    end
end